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ABSTRACT 

We present a comprehensive spectral and morphological analysis of the remnant of Supernova (SN) 
1987A with the Australia Telescope Compact Array (ATCA) and the Atacama Large Millime- 
ter/submillimeter Array (ALMA). The non-thermal and thermal components of the radio emission 
are investigated in images from 94 to 672 GHz (A 3.2 mm to 450 /im), with the assistance of a 
high-resolution 44 GHz synchrotron template from the ATCA, and a dust template from ALMA ob- 
servations at 672 GHz. An analysis of the emission distribution over the equatorial ring in images 
from 44 to 345 GHz highlights a gradual decrease of the east-to-west asymmetry ratio with frequency. 
We attribute this to the shorter synchrotron lifetime at high frequencies. Across the transition from 
radio to far infrared, both the synchrotron/dust-subtracted images and the spectral energy distribu- 
tion (SED) suggest additional emission beside the main synchrotron component (S u ex v~ 0 - 73 ) and the 
thermal component originating from dust grains at T ~ 22 K. This excess could be due to free- free 
flux or emission from grains of colder dust. However, a second flat-spectrum synchrotron component 
appears to better fit the SED, implying that the emission could be attributed to a pulsar wind nebula 
(PWN). The residual emission is mainly localised west of the SN site, as the spectral analysis yields 
—0.4 < a < —0.1 across the western regions, with a ~ 0 around the central region. If there is a 
PWN in the remnant interior, these data suggest that the pulsar may be offset westward from the SN 
position. 

Keywords: radio continuum: general — supernovae: individual (SN 1987A) — ISM: supernova 
remnants — radiation mechanisms: non-thermal — radiation mechanisms: thermal — 
stars: neutron 



1. INTRODUCTION 

The evolution of the remnant of Supernova (SN) 1987A 
in the Large Magellanic Cloud has been closely mon- 
itored since the collapse of its progenitor star, Sand- 
uleak (Sk) -69° 202, on 1987 February 23. Models of 
Sk —69° 202 indicated that it had an initial mass of ~20 



Mq (Hillebrandt et al. 1987). The mass range of the 
progenitor is consistent with the formation of a neutron 
star (Thielemann & Arnett 1985), and thus with the 
neutrino events reported by the KamiokaNDE (Hirata 
et al. 1987) and 1MB (Bionta et al. 1987; Haines et 
al. 1988) detectors. Models by Crotts & Heathcote 
(2000) suggest a transition from red supergiant (RSG) 
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into blue supergiant (BSG) to explain the hourglass neb- 
ula structure, which envelopes the SN with three nearly- 
stationary rings (Blondin & Lundqvist, 1993; Chevalier 
& Dwarkadas 1995; Martin & Arnett 1995; Morris & 
Podsiadlowski 2007). The two outer rings, imaged with 
the Hubble Space Telescope (HST; Jakobsen et al. 1991; 
Plait et al. 1995), are located on either side of the cen- 
tral ring in the equatorial plane (equatorial ring, ER), 
and likely formed at the same time as the ER (Crotts 
& Heathcote 2000; Tziamtzis et al. 2011). The syn- 
chrotron emission, generated by the shock propagating 
into the clumpy circumstellar medium (CSM) close to the 
equatorial plane, was detected in the mid-90s (Staveley- 
Smith et al. 1992; Turtle et al. 1990), and has become 
brighter over time (Manchester et al. 2002; Zanardo et 
al. 2010). Radio observations have stretched from flux 
monitoring at 843 MHz with the Molonglo Observatory 
Synthesis Telescope (Ball et al. 2001; Staveley-Smith 
et al. 1993) to images of sub-arcsec resolution with the 
Australia Telescope Compact Array (ATCA) (Gaensler 
et al. 1997; Manchester et al. 2005; Ng et al. 2008; 
Potter et al. 2009; Zanardo et al. 2013). ATCA ob- 
servations at 94 GHz (Lakicevic et al. 2012a) have been 
followed by observations from 100 GHz up to 680 GHz 
with the Atacama Large Millimeter /submillimeter Array 
(ALMA; Indebetouw et al. 2014; Kamenetzky et al. 
2013). 

The ongoing shock expansion has been monitored at 9 
GHz since 1992 (Gaensler et al. 1997, 2007). Shock 
velocities of ~4000 km s _1 have been measured between 
day 4000 and 7000 (Ng et al. 2008), while signs of 
a slower expansion have been tentatively detected after 
day ~7000, as the shock has likely propagated past the 
high-density CSM in the ER (Ng et al. 2013). Simi- 
lar evidence of slower shock expansion since day ~6000 
has been found in X-ray data (Park et al. 2005, 2006; 
Racusin et al. 2009) as well as in infrared (IR) data 
(Bouchet et al. 2006). 

Since the early super-resolved images at 9 GHz (Gaensler 
et al. 1997), a limb-brightened shell morphology has 
been characteristic of the remnant. The radio emission, 
over the years, has become more similar to an elliptical 
ring rather than the original truncated-shell torus (Ng 
et al. 2013). The radio remnant has shown a consis- 
tent east- west asymmetry peaking on the eastern lobe, 
which has been associated with higher expansion veloc- 
ities of the eastbound shocks (Zanardo et al. 2013). 
The asymmetry degree appears to have changed with the 
shock expansion, as images at 9 GHz exhibit a decreasing 
trend in the east-west asymmetry since day ^7000 (Ng et 
al. 2013). High-resolution observations at 1.4-1.6 GHz 
(Ng et al. 2011) via Very Long Baseline Interferometry 
(VLBI) with the Australian Large Baseline Array (LB A), 
have highlighted the presence of small-scale structures in 
the brightest regions in both lobes. 

The relation between the radio emission and the syn- 
chrotron spectral indices, a (S u oc has been inves- 
tigated via both flux monitoring (Manchester et al. 
2002; Zanardo et al. 2010) and imaging observations 
(Lakicevic et al. 2012a; Potter et al. 2009; Zanardo 
et al. 2013) with the ATCA. The progressive flatten- 
ing of the radio spectrum derived from 843 MHz to 8.6 



GHz at least since day 5000, coupled with the e-folding 
rate of the radio emission, has pointed to an increas- 
ing production of non-thermal electrons and cosmic rays 
(CR) by the shock front (Zanardo et al. 2010). On the 
other hand, the association of steeper spectral indices 
with the brightest eastern sites implies a higher injection 
efficiency on the eastern side of the SNR (Zanardo et al. 
2013). Flatter spectral indices in the center of the rem- 
nant have been tentatively identified from low-resolution 
two-frequency spectral maps (Lakicevic et al. 2012a; 
Potter et al. 2009), while at 18-44 GHz the central and 
center-north regions have —0.5 < a < —0.3 (Zanardo et 
al. 2013). With ALMA, the spectral energy distribution 
(SED) of the remnant has been mapped where the non- 
thermal and thermal components of the emission overlap, 
identifying cold dust in the SNR interior (Indebetouw et 
al. 2014), which accommodates 0.4 — 0.7 M 0 of the dust 
mass discovered with the Herschel Space Observatory 
(Herschel; Matsuura et al. 2011) in the ejecta. 

This paper combines the results presented by Indebetouw 
et al. (2014) with a comprehensive morphological and 
spectral analysis of SNR 1987A based on both ATCA 
and ALMA data. In § 2, we present the ALMA Cycle 

0 super-resolved images before and after subtraction, in 
the Fourier plane, of the synchrotron and dust compo- 
nents (§ 3). In § 4, we assess the remnant asymmetry 
from 44 to 345 GHz. In § 5, we update the SED derived 
by Indebetouw et al. (2014), while, in § 6, we investigate 
the spectral index variations in the SNR across the tran- 
sition from radio to far infrared (FIR). In § 7, we discuss 
possible particle flux injection by a pulsar situated in the 
inner regions of the remnant. 

2. OBSERVATIONS AND ANALYSIS 

The ATCA and ALMA observations used in this study 
were performed in 2011 and 2012. ATCA observations 
at 44 and 94 GHz are detailed in Zanardo et al. (2013) 
and Lakicevic et al. (2012a), respectively. ALMA ob- 
servations were made in 2012 (Cycle 0) from April to 
November, over four frequency bands: Band 3 (B3, 84- 
116 GHz, A 3 mm), Band 6 (B6, 211-275 GHz, A 1.3 
mm), Band 7 (B7, 275-373 GHz, A 850 /im) and Band 
9 (B9, 602-720 GHz, A 450 /im). Each band was split 
over dual 2-GHz-wide sidebands, with minimum base- 
lines of 17 m (B9) to maximum baselines of 400 m (B3). 
All observations used quasars J0538-440 and J0637-752 
as bandpass and phase calibrators, respectively. Callisto 
was observed as an absolute flux calibrator in B3 and B6, 
while Ceres was used in B7 and B9 (see also Kamenetzky 
et al. 2013). It is noted that, while ALMA is designed 
to yield data with flux density calibration uncertainty as 
low as ~1%, in Cycle 0 this uncertainty is estimated at 
~5% at all frequencies. Relevant observational parame- 
ters are listed in Table 1 (see also Table 1 in Indebetouw 
et al. 2014). 

Each dataset was calibrated with the CAS A 1 package, 
then exported in miriad 2 for imaging. After clean- 
ing (Hogbom 1974), both phase and amplitude self- 
calibration were applied in B3 over a 2-minute solution 

1 http: / /casa. nrao.edu / 

2 http: / / www. atnf .csiro. au / computing/software / miriad / 
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Table 1 

ALMA Observing Parameters 



Parameter 


i no Ur, 


Zio LjrrlZ 


O40 <or IT Z 


Of Z vjtITZ 




(Band 3, B3) 


(Band 6, B6) 


(Band 7, B7) 


(Band 9, 55) 


Date (2012) 


Apr 5, 6 


Jul 15 & Aug 10 


Jul 14 & Aug 24 


Aug 25 & Nov 5 


Day since explosion 


9174 


9287 


9294 


9351 


Frequency bands* (GHz) 


100.093-101.949 


213.506-213.597 


336.979-340.917 


661.992-665.992 


102.051-103.907 




349.010-352.963 


678.008-682.008 


Center frequency (GHz) 


101.918 


213.146 


345.364 


672.165 


Channel width (MHz) 


4.883 


4.883 


31.250 


15.625 


Max baselines (u, v) (kA) 


150, 120 


180, 260 


400, 400 


700, 700 


No. of antennas 


14-18 


14-23 


28 


19-25 


Total observing time (hr) 


0.83 


1.03 


0.62 


3.40 



* In B3 and B6, the frequency range is selected to avoid CO and SiO emission (Kamenetzky et al. 
2013). 



Table 2 

Image Parameters 



Image 


s v (-) 


SR Beam ( b ) 


DL Beam ( c ) 


PA 


Rms noise 


Dynamic range 


(GHz) 


(mJy) 


(") 


(") 


(°) 


(mjy/beam) 




94 (d) 


24.2 ±3.9 


0.7 






0.085 


137 


102 (e) 


23.1 ±3.1 


0.8 


1.74 x 1.25 


5.7 


0.033 


285 


213 


19.7 ± 1.6 


0.6 


1.16 x 0.74 


-68.5 


0.034 


75 


345 


16.7 ± 1.5 


0.3 


0.65 x 0.48 


-40.2 


0.023 


121 


672 


52.8 ± 14.2 


0.3 


0.34 x 0.28 


68.5 


1.219 


28 


94-/ B 9 (/) 


23.2 ±3.9 


0.7 


0.78 x 0.63 


15.4 


0.094 


81 


102-/B9 


19.4 ±3.2 


0.8 


1.66 x 1.19 


6.3 


0.034 


143 


213-/B9 


16.9 ±1.9 


0.6 


1.16 x 0.74 


-68.5 


0.029 


64 


345-/B9 


11.5 ±1.9 


0.3 


0.65 x 0.48 


-44.0 


0.031 


43 


672-/B9 


-1.0 ±14.3 


0.3 


0.34 x 0.28 


68.5 


1.193 


2 


94-/ 4 4 (9) 


0.9 ±3.9 


0.7 


0.78 x 0.63 


15.4 


0.095 


60 


102-/44 


3.5 ±3.1 


0.8 


1.48 x 0.96 


10.1 


0.027 


81 


213-/44 


2.9 ± 1.7 


0.6 


1.16 x 0.74 


-68.5 


0.024 


26 


345-/44 


5.8 ± 1.6 


0.3 


0.65 x 0.48 


-44.0 


0.014 


90 


672-/44 


47.4 ± 14.8 


0.3 


0.57 x 0.55 


78.6 


0.816 


24 



a All flux densities are derived from the diffraction-limited images. The errors are derived from the flux 

calibration uncertainty combined with the uncertainties in the image subtraction. 

b Circular beam used for super-resolution (SR). 

c Beam associated with the diffraction-limited (DL) image. 

d The flux density is scaled to day 9280 via exponential fitting parameters derived for ATCA flux densities 
from day 8000, as measured at 8.6 and 9 GHz (Zanardo et al., in preparation). All other image parameters 
are as from Lakicevic et al. (2012a). 
e Images are shown in Figures 1 and 2. 

f Images obtained by subtracting the model flux density at 672 GHz (/eg) scaled to fit the central emission. 
See central column of Figures 1 and 2. 

9 Images obtained by subtracting the model flux density at 44 GHz (I 44) scaled to fit the toroidal emission. 
See left column of Figures 1 and 2. 



interval, while only phase calibration was applied in B6 
and B7. No self-calibration was performed in B9. As 
in Zanardo et al. (2013), we note that since the self- 
calibration technique removes position information, each 
image was compared with that prior to self-calibration 
and, in case of positional changes, the self-calibrated im- 
ages were shifted. Further adjustments were made in 
the comparison with the ATCA observations at 44 GHz, 
based on prominent features on the eastern lobe and lo- 
cation of the remnant center. As from Zanardo et al. 
(2013), the 44 GHz image was aligned with VLBI ob- 
servations of the SNR (Zanardo et al. in preparation). 
Adding in quadrature these positional uncertainties and 
the accuracy of the LB A VLBI frame, the errors in the 
final image position are estimated at ~60 mas. 



Deconvolution was carried out via the maximum entropy 
method (MEM) (Gull & Daniell 1978) in B3, B6 and B7. 
A weighting parameter of robust = 0.5 (Briggs 1995) 
was used in all bands. The resultant diffraction-limited 
images, which have central frequency at 102 GHz in B3, 
213 GHz in B6, 345 GHz in B7 and 672 GHz in B9, 
were then super-resolved with a circular beam of 0'/8 in 
B3, Of! 6 in B6, and 0'/3 in B7 and B9. The diffraction 
limited and super-resolved images are shown in the first 
column of Figures 1 and 2, below the ATCA image at 
94 GHz (Lakicevic et al. 2012a). Integrated Stokes / 
flux densities, dynamic range and related rms are given 
in Table 2. 

To decouple the non-thermal emission from that orig- 
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Figure 1. Top to bottom - Left column: Stokes / continuum images of SNR 1987A at 94 (Lakicevic et al. 2012a), 102, 213, 345 and 
672 GHz. Images from 102 to 672 GHz are made from ALMA observations (Cycle 0) performed from 2012 April 5 to November 5 (see 
Table 1). center column: Images obtained by subtracting the model flux density at 672 GHz (Band 9) scaled to fit the central emission. 
Right column: Images obtained by subtracting a scaled model flux density at 44 GHz (Zanardo et al. 2013), with 3<r flux density contours 
highlighted (white). The angular resolution is shown in the bottom left corner. The green cross indicates the VLBI position of SN 1987A 
as determined by Reynolds et al. (1995) [RA 05 h 35 m 27!968, Dec -69° 16' 11//09 (J2000)]. 
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28.20s 28.00s 5h35m27.80s 

RA (J 2 000) 



28.20s 28.00s 5h35m27.80s 

RA (J2000) 



28.20s 28.00s 5h35m27.80s 

RA (J 2 000) 



Figure 2. The panel layout is identical to that in Figure 1, but all Stokes I continuum subtracted images are super-resolved. The circular 
beam used for the super-resolved ALMA images is 0'/8 at 102 GHz (Band 3), 0'/6 at 213 GHz (Band 6), Of.' 3 at 345 (Band 7) and 672 
GHz (Band 9), and is plotted in the lower left corner. The 94 GHz images are restored with a circular beam. For images on the right 
column, 3cr flux density contours are highlighted (white). The green cross indicates the VLBI position of SN 1987A as determined by 
Reynolds et al. (1995). 
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Figure 3. Stokes I continuum image of SNR 1987A at 213 GHz as obtained after the dual subtraction of the model flux densities 
at 44 and 672 GHz (see § 2). The image is super-resolved with a CK'6 circular beam, and the yellow contour highlights the 3cr flux 
density levels. The off-source rms noise is 0.14 mJy beam -1 (S/N ~ 6 a). For comparison, the 213 GHz image obtained after the single 
subtraction of the model flux density at 44 GHz and similarly super-resolved with a Of/ 6 circular beam (see Figure 2), is outlined via the 
3cr flux density levels (magenta contours). The model images, i.e. the Stokes I images at 44 GHz (Zanardo et al. 2013) and at 672 GHz, 
are outlined by blue and white contours, respectively, at the 5 a flux density levels. Both model images are resolved with a Of/ 3 circular beam. 



inating from dust, the synchrotron component, as re- 
solved with ATCA at 44 GHz (Zanardo et al. 2013), 
and the dust component, as imaged with ALMA at 672 
GHz (B9) (Indebetouw et al. 2014), were separately sub- 
tracted from the datasets at 94, 102, 213, 345 and 672 
GHz. All subtractions were performed in the Fourier 
plane, via MIRIAD task uvmodel, where the model flux 
density at 44 GHz was scaled to fit the SNR emission 
over the ER (I44), while the B9 model flux density was 
scaled to fit the emission localised in the central region 
of the remnant (2b9)- Scaling of the 44 GHz model was 
tuned by minimizing the flux density difference on the 
brighter eastern lobe, without over-subtracting in other 
regions of the remnant. To separate the emission in the 
SNR center, the central flux was firstly estimated by fit- 
ting a gaussian model via miriad task uvf it. The image 
model at 672 GHz was then scaled to match the flux of 
the gaussian model. The scaling factor was further tuned 
to minimize over-subtraction. 

Super-Nyquist sampling was applied in all images, using 
a pixel size of 8 mas to avoid artefacts when sources are 
not at pixel centers. Deconvolution via MEM was car- 
ried out on the residual images obtained from the sub- 
traction of Jb9, while standard CLEANing was applied 
to the residuals obtained from the I 44 subtraction. All 
diffraction-limited subtracted images are shown in the 
central and right columns of Figure 1, while Figure 2 
shows the residuals after super-resolution with the circu- 
lar beam used for the original images. All image param- 



eters are given in Table 2. 

The flux densities were determined by integrating within 
polygons enclosing the SNR emission. Uncertainties in 
the flux densities include uncertainties in the image fit- 
ting/scaling process combined with the uncertainty in 
the flux density calibration. We note that the residual 
images from the subtraction of both models at 44 and 
672 GHz were not considered, since the error attached to 
the double subtraction exceeds the total integrated flux 
density. Only at 213 GHz the error-flux margin is minor, 
as the diffraction-limited image obtained after the dual 
subtraction has integrated flux density of 2.1 ± 1.9 mJy 
(Figure 3). 

3. MORPHOLOGY 

ALMA observations of SNR 1987A capture both the rem- 
nant emission from the ER and that from the SNR inte- 
rior, where the dense ejecta sit (Indebetouw et al. 2014; 
Kamenetzky et al. 2013). 

While the image at 102 GHz barely resolves the two- 
lobe distribution of the emission, the images at 213 and 
345 GHz clearly show the ringlike emission morphology, 
localised around the ER (see Figures 1 and 2). It is un- 
derstood that the radio emission over the ER is primarily 
synchrotron emission, generated by the interaction of the 
SN shock with the dense CSM near the equatorial plane, 
which results in a magnetic-field discontinuity where par- 
ticles are accelerated (e.g. Zanardo et al. 2010, 2013). 
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Figure 4. Overlay of the HST image of SNR 1987A (light blue contours; (Larsson et al. 2011)) on the 345 GHz image produced from 
ALMA observations performed in 2012 June and August (brown-yellow color scale and red contours). As from Figure 2, the ALMA image 
is super-resolved with a (y/3 circular beam. 



As shown in early models of Type II SNRs (Chevalier 
1982), the region of interaction between the SN blast 
and the CSM consists of a double-shock structure, with 
a forward shock, where ambient gas is compressed and 
heated, and a reverse shock, where the ejecta are decel- 
erated. Between the two shocks, the reflected shocks, 
due to the forward shock colliding with the dense ER 
(Borkowski, Blondin & McCray 1997), propagate inward 
(Zhekov et al. 2009, 2010). The reverse shock at first 
expanded outwards, behind the forward shock, but might 
have been inward- moving since day ~ 7000. As discussed 
by Ng et al. (2013), the ring synchrotron emission is cur- 
rently localised between the forward and reverse shocks, 
and likely has components from both the ER and high- 
latitude material above the equatorial plane. Truncated- 
shell torus models of the remnant geometry at 9 GHz 
indicate that the half-opening angle has been decreasing 
since day ~7000, and is estimated at ~ 27° at day 9568 
(Ng et al. 2013). 

At 345 GHz the SNR interior is brighter, while at 672 
GHz the emission is predominantly localised in the cen- 
tral region of the remnant. Since the emission from this 
region rises steeply with frequency, it has been identified 
with thermal dust emission (Indebetouw et al. 2014), 
as dust grains, probably heated by 44 Ti decay and X-ray 
emission from the reverse shock (Larsson et al. 2013), 
emit strongly in the FIR regime. The central emission, 
visible both in B7 and B9, appears to extend over the in- 
ner optical ejecta (see Figure 4). In particular, as noted 



by Indebetouw et al. (2014), this inner emission shows a 
north-south elongation, which in B9 can be identified be- 
tween PA ~ 20° and PA ~ 30°, similar to that seen with 
HST (Larsson et al. 2011, 2013). From both Figures 1 
and 2, it can also be noticed that the SNR emission at 
672 GHz includes possible emission located to the NW 
(see contour overlays in Figure 5). This NW feature has 
signal-to-noise ratio (S/N) of 3.7 a, and integrated flux 
density of 6.8d=0.5 mJy, i.e. ~10% of the total integrated 
flux density at 672 GHz. 

3.1. Subtracted Images 

To identify the origin of the emission in the 
244— subtracted images, with respect to the structure of 
the remnant as seen with HST (Larsson et al. 2011), the 
diffraction-limited residuals at 102, 213 and 345 GHz are 
superimposed in Figure 5. It can be seen that the resid- 
ual at 102 GHz, characterised by S/N| 102 = 6.9 <r, is 
mainly located on the western lobe, west of the VLBI 
position of SN 1987A, as determined by Reynolds et 
al. (1995) [RA 05 h 35 m 27!968, Dec -69° 16' IV! 09 
(J2000)]. We note that, since the residual emission at 
102 GHz has a constant flux per synthesized beam, the 
S/N ratio is independent on the beam size. The emis- 
sion at 213 GHz (orange contours in Figure 5), with 
S/N| 213 = 7.7 a, peaks NW of the SN position, while 

fainter emission (S/N| 2 5 13 = 4.2 a) may extend NE. The 
residual emission at both 102 and 213 GHz is above 
noise levels and, thus, unlikely to be the result of image 
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Figure 5. Comparison of the diffraction-limited images of SNR 1987A at 102, 213 and 345 GHz, as obtained after subtraction of the 
scaled model flux density at 44 GHz (see Figure 1, right column), with the optical image of the remnant. In detail, the residual image 
at 102 GHz (brown-yellow color scale) is overlaid with the contours outlining the residual images at 213 GHz (orange) and at 345 GHz 
(yellow). To locate the main sites associated with dust emission, the diffraction-limited image at 672 GHz is also outlined (red) (as from 
Figure 1, left column). Images at 213, 345 and 672 GHz are overlaid with contours at 3, 4 and 5<r flux density levels. The angular beams 
associated with each image, as shown in Figure 1, are listed in Table 2. The HST image (Larsson et al. 2011) is outlined via contours 
(blue) that highlight the structure of the outer rings, the equatorial ring and the ejecta. 



artefacts. Given the brightening of the emission from 
the dust in the central region of the SNR, the resid- 
ual images at 345 and 672 GHz have higher S/N. In 
particular, at 345 GHz (yellow contours in Figure 5) 
S/N| 345 = 17.4cr, while, similar to the morphology at 
213 GHz, the residual emission extends westwards and 
elongates NW, with a much fainter spot on the north- 
eastern section of the ER ( S/N|^ 5 = 2.4 a). A westward- 
elongated morphology is present in the image at 672 GHz 
(S/N| 672 = 12.5 a). 

The Ib9 _ subtracted images emphasize the ringlike mor- 
phology of the synchrotron emission that mainly origi- 
nates near the ER (see § 3) and, thus, the asymmetry of 
this emission, which is discussed in § 4. These residuals 
also highlight the presence of NW emission at 672 GHz, 
i.e. outside the inner SNR (see red contours in Figure 5). 
In fact, at 94, 213, and 345 GHz, it can be seen that 



the discontinuity in the NW sector of the ER, between 
PA ~ 290° and PA ~ 300°, becomes more prominent 
after subtraction of the B9 model. Some extended emis- 
sion north and NW of the ER emerges from noise at 345 
GHz. 

4. ASYMMETRY 

The east-west asymmetry of the synchrotron emission, 
primarily associated with the emission morphology over 
the ER, is investigated from 44 to 345 GHz, as shown 
in Figure 6, where all images are restored with a (X'7 
circular beam. At 213 and 345 GHz the central dust 
emission is subtracted, via a scaled model of the flux 
density as resolved in B9 (see § 2). 

Considering radial sections crossing the SN site 
(Reynolds et al. 1995), the SNR asymmetry is firstly 
estimated as the ratio between the eastern peak of the 
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Figure 6. Top to bottom - Left column: Stokes I continuum images of SNR 1987A at 44, 94, 102, 213, and 345 GHz. The 213 (213-B9) 
and 345 GHz (345-B9) images are derived after subtraction of the scaled model flux density at 672 GHz (Band 9, B9). All images are 
resolved with a 0 x /7 circular beam (as shown in the bottom left corner of the 44 GHz image). Center column: The images on the left 
column are converted to polar coordinates to visualize asymmetries in the two-dimensional radial distribution. The position angle and the 
projected radial distance from the geometrical center of the remnant, are the new coordinates. Each image is projected along concentric 
ellipses (e.g. blue ellipse overlaid on the 44 GHz image). The emission intensity is conserved in the conversion. Contours at 65% and 85% 
levels of the emission are shown on all images {white). Right column: Radial slices through each image at four position angles (as indicated 
on the images on the left column), which correspond to 0° {light green), 90° {dark green), emission peaks on the eastern lobe {white/black) 
and on the western lobe {red). The radial slices are plotted against the RA offset from the VLBI position of SN 1987A as determined by 
Reynolds et al. (1995) [RA 05 h 35 m 27!968, Dec -69° 16' ll 7 / 09 (J2000)]. 
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radial slice crossing the maximum emission on the east- 
ern lobe (see black/white profile in Figure 6), and the 
western peak of the radial slice crossing the maximum 
emission on the western lobe, (see red profiles in Fig- 
ure 6). This ratio, A p , emphasizes the hot spots on each 
side of the remnant. Alternatively, the ratio Ai between 
the total flux densities integrated over the eastern and 
western halves of the image, is derived by splitting the 
image at the RA associated with the SN site. For com- 
parison with previous asymmetry estimates (e.g. Zanardo 
et al. 2013), the ratio A§§ between the eastern and west- 
ern peaks of the radial slice at PA 90° (see dark green 
profile in Figure 6) is also derived. All values for A p , A{ 
and Ago, both for the images restored with a 0'/7 circular 
beam (as in Figure 6) and the diffraction-limited images, 
are plotted in Figure 7 and listed in Table 3. 

From Figure 6 it can be seen that the Ago ratio does 
not fully capture the asymmetry changes with frequency, 
since the hot spot in the western lobe, which becomes 
brighter from 102 to 345 GHz, is located southwards of 
the 90° profile. In Figure 7, the linear fits derived for 
both A p and Ai ratios show a consistent decrease as fre- 
quencies reach the FIR. At 345 GHz, Ai values indicate 
that the east- west asymmetry is reversed, thus matching 
the asymmetry trend seen in recent HST images (Lars- 
son et al. 2011), where the western side of the ring is 
markedly brighter. As discussed in § 3, the morphology 
similarities between the optical image of the SNR and 
the super-resolved image at 345 GHz are evident (see 
Figure 4). 

The change of the remnant's east- west asymmetry over 
time has been discussed by Ng et al. (2013), as the 
result of a progressive flattening of the shock structure 
in the equatorial plane, due to the shock becoming en- 
gulfed in the dense UV-optical knots in the ER, coupled 
with faster shocks in the east side of the remnant. While 
X-ray observations do not show significant difference be- 
tween the NE and SW reverse shock velocities, although 
the NE sector is brighter (Fransson et al. 2013), faster 
eastbound outer shocks have been measured in the radio 
(Zanardo et al. 2013) and point to an asymmetric ex- 
plosion of a binary merger as SN progenitor (Morris & 
Podsiadlowski 2007, 2009). As the SN blast is gradually 
overtaking the ER, faster expanding shocks in the east 
would exit the ER earlier than in the west. 

The effects of the asymmetric shock propagation are 
likely to emerge in the transition from radio to FIR 
rather than at lower frequencies, due to the shorter syn- 
chrotron lifetime at higher frequencies. To estimate the 
synchrotron lifetime in the FIR range, we use the ap- 
proximation that, in a magnetic field of strength B, all 
the radiation of an electron of energy E is emitted only 
at the critical frequency v c (Rybicki & Lightman 1979). 
Considering the electrons orbit is inclined at a pitch angle 
0 to the magnetic field, the synchrotron lifetime, r e , can 
be derived as a function of v c (e.g. Condon 1992) 

Te = WJdt\ ~ L06 x lo9 ( Bsin #)" 3/ V 1/2 , (i) 

where r e is expressed in years, B in fiQ and v in GHz. 
As the electrons are expected to have an isotropic distri- 
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Figure 7. Asymmetry ratios of Stokes I images of SNR 1987A 
from 44 to 345 GHz. With reference to Table 3, the east-west 
asymmetry ratios estimated via three different methods are 
indicated as: A p (black symbols); Ai (red symbols); Ago (blue 
symbols). Asymmetry ratios are derived from the diffraction- 
limited images (diamonds) (see Figure 1), from the images 
restored with a u/7 circular beam (circles) (see Figure 6); and 
from images obtained after the subtraction of the dust emission 
component (stars) (see Figure 6). Tentative linear fits are shown 
for A p (black fits), Ai (red fits) and Ago (blue fits), based on 
all ratios derived via each method. The orthogonal distance 
regression (Aodr) and the robust fitting from the square root 
of absolute residuals (A^bs), are tested for linear interpolation. 
The combined fits for the A p and A{ ratios yield an average 
A(v) nu A{vq) — 1.9 db 0.8 x 10~ 3 v, where v is expressed in GHz 
and A(vq = 30) = 1.57 ± 0.13. 

bution of pitch angles, (sin 2 (9) = 2/3; for 200 < v < 400 
GHz, and assuming a magnetic field strength at the 
shock front of 10 < B < 20 mG (Berezhko & Kseno- 
fontov 2006; Berezhko et al. 2011), in the radio/FIR 
transition we can estimate 20 < r e < 80 yr. However, 
since it is likely that sub-diffusive particle transport 
(Kirk et al. 1996) is taking place in regions of the SNR, 
in conjunction with efficient CR acceleration (Glushkov 
et al. 2013), local magnetic-field amplifications could 
exceed the above limits by at least an order of mag- 
nitude (Bell & Lucek 2001). In this scenario, the 
local r e would be of the order of months. Therefore, 
with regions in the remnant where electrons might be 
unable to cross the emission sites within their radiative 
life- time, the synchrotron emission at FIR frequencies 
would require the presence of relatively fresh injected 
and/or re- accelerated electrons to match the emission 
distribution at lower frequencies. 

5. SPECTRAL ENERGY DISTRIBUTION 

The spectral energy distribution for ATCA data from 
I A to 94 GHz and ALMA data is shown in Figure 8. To 
match the average epoch of the ALMA data, the ATCA 
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Table 3 

Asymmetry Ratios 



y A _ S u (E Max ) (a) _ r / r (b) —of E Max \ 

Image A p = — — A t = j S u / j S v A 90 = S u — 



GHz 




0"7 (d) 


DL (e) 


0"7 


DL 


0'/7 


DL 


44 




1.61 ±0.08 


1.42 ±0.05 


1.36 ±0.02 


1.53 ±0.07 


1.69 ±0.08 


1.54 ±0.05 


94 




1.55 ±0.05 




1.18 ±0.04 




1.88 ±0.07 




102 




1.44 ±0.04 




1.10 ±0.04 




1.40 ±0.07 




213 




1.28 ±0.03 


1.43 ±0.05 


1.24 ±0.06 


1.23 ±0.02 


1.43 ±0.02 


1.49 ±0.05 


213- 


W /} 


1.44 ±0.03 




1.30 ±0.02 




1.19 ±0.03 




345 






1.28 ±0.03 




0.78 ±0.04 




1.16 ±0.03 


345- 


7 B 9 (/) 


1.30 ±0.02 




0.80 ±0.02 




1.48 ±0.04 





a Ratio between the eastern peak of the radial slice crossing the maximum emission on the eastern 
lobe, and the western peak of the radial slice crossing the maximum emission on the western lobe 
(see black/white and red profiles in Figure 6). 

b Ratio between the total flux densities integrated over the eastern and western halves of the image. 
The image is split at the RA associated with the SN site (Reynolds et al. 1995). 
c Ratio between the eastern and western peaks of the radial slice at PA 90° (see dark green profile 
in Figure 6). 

d Images resolved with a (Y/7 circular beam (see Figure 6). 
e Diffraction-limited (DL) images (see Figure 1). 

f Images derived after subtraction of the model flux density at 672 GHz (Ib9), scaled to fit the 
central emission (see Figure 6). 



data are scaled to day 9280, via exponential fitting pa- 
rameters derived for ATCA flux densities from day 8000, 
as measured at 1.4, 8.6 and 9 GHz (Zanardo et al., in 
preparation; Staveley-Smith et al. 2014). Across the 
transition from radio to FIR frequencies, the observed 
spectrum consists of the sum of thermal and non-thermal 
components. 

As described in § 2, to identify the dust component of the 
emission from the inner regions of the remnant, the B9 
model flux density, 5b9, has been scaled to fit the emis- 
sion measured in the SNR central region at 94-345 GHz 
(*SB9fit> hollow red circles in Figure 8). The subtraction 
of SB9 fit from the visibilities at 94-345 GHz yields the 
residual flux densities indicated as S v — $B9 fit (purple 
bars in Figure 8), as for the images shown in the central 
column of Figure 1 . The flux densities Sb9 fit derived in 
B6 and B7, together with the total integrated flux density 
measured in B9 (Sb9), although obtained via a different 
reduction technique, have been associated by Indebetouw 
et al. (2014) with dust grains, in conjunction with data 
from Herschel (Matsuura et al. 2011) and the Atacama 
Pathfinder Experiment (APEX; Lakicevic et al. 2012b) 
(see Figure 8). 

Similarly (see § 2), to separate the non-thermal emis- 
sion from that thermal, the 44 GHz model flux den- 
sity, S44, has been scaled to fit the toroidal compo- 
nent of the emission at 94-672 GHz (blue/cyan dia- 
monds in Figure 8). By fitting the resulting 644 fit com- 
ponents and the ATCA flux densities at 1.4-44 GHz, 
we obtain, at day 9280, the synchrotron spectral index 
a s = -0.727 ± 0.020, with S u (a s ) oc v as . The spectral 
index measured from 1.4 to 94 GHz, i.e. for ATCA data 
only, is a ATCA = -0.735±0.028. While a s is slightly flat- 
ter than a ATCA , both values are consistent with the pro- 



gressive flattening of the radio spectrum measured since 
day ~ 5000 (Zanardo et al. 2010). The subtraction 
of 644 fit from the visibilities at 94-672 GHz, gives the 
residual flux densities S v — S±± fit (red/orange diamonds 
in Figure 8), as for the images shown in the right column 
of Figures 1. In B3, B6 and B7, S v - 5 4 4 fit > 5B9 fit , 
i.e. the residuals exceed the emission expected from the 
dust. 

While the subtraction of the flux densities is in- 
evitably affected by errors (see Table 2), given that the 
244— subtracted images have a S/N> 7 a and the resid- 
ual emission appears primarily located westwards of the 
optical ejecta (see § 3.1), we investigate the nature of 
this emission excess as: (1) free- free emission from an 
ionized fraction of the inner ejecta; (2) synchrotron emis- 
sion from a compact source located in the inner regions 
of the remnant; (3) emission from grains of very cold 
dust. 

5.1. Free-Free Emission 

To estimate the free-free radiation in the SNR as im- 
aged with ALMA, we hypothesize an ionized portion 
of the ejecta as an approximately spherical region, lo- 
cated inside the ER. Using the beam size of the super- 
resolved ^44— subtracted image at 102 GHz as an upper 
limit, we consider the radius of the spherical region up 
to R 8 - 0'/40 (w 3.08 x 10 17 cm). Such radius covers 
the extent of the inner ejecta as imaged in the optical 
(Larsson et al. 2013), and stretches to the likely ra- 
dius of the reverse shock, qualitatively identified with 
the inner edge of the emission over the ER (see Fig- 
ure 4). Given that the pre-supernova mass has been 
estimated between 14 M 0 and 20 M 0 (Smartt 2009), 
if we assume that the Hll region in the ejecta has a uni- 
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Figure 8. Spectral energy distribution (SED) of SNR 1987A from radio to FIR, with data from: ATCA at 1.4 GHz (Zanardo et al. in 
prep), 9 GHz (Ng et al. 2013), 18 and 44 GHz (Zanardo et al. 2013), and 94 GHz (Lakicevic et al. 2012a); ALMA at 102 GHz (Band 
3, B3), 213 GHz (Band 6, B6), 345 GHz (Band 7, B7) and 672 GHz (Band 9, B9) (this paper); the Atacama Pathfinder Experiment 
(APEX) at 345 and 857 GHz (Lakicevic et al. 2012b); the Herschel Space Observatory at 600 - 3000 GHz (Matsuura et al. 2011). The 
brown dash-dot-dotted curve is the amorphous carbon dust fit for ALMA data and Herschel observations carried out in 2012 (Matsuura 
et al. 2014, submitted). To match the average epoch of ALMA observations, ATCA data are scaled to day 9280, via exponential fitting 
parameters derived for ATCA flux densities from day 8000, as measured at 1.4, 8.6 and 9 GHz (Zanardo et al., in preparation; Staveley- 
Smith et al. 2014). The hollow red diamond indicates the central emission measured at 44 GHz as reported by Zanardo et al. (2013), 
scaled to day 9280. The difference between possible spectral fits is highlighted in light orange. See § 5 for a detailed description of this figure. 



form density p e j( r ) = 3MHii(47rr 3 ) x , we set 0.7M© < 
Mhii < 2.5 M©. The lower limit, M H n ~ 0.7 M©, i.e. 
p ej (r = R s ) = 1.1 x 10 -20 , represents a partial ionization 
of the ejecta by X-ray flux, either within the inner region 
or on the outer layer. The upper limit, Mhii ~ 2.5 M©, 
i.e. p e j(r = R s ) = 4.1 x 1CT 20 g cm -3 matches the 
density model by Blinnikov et al. (2000) scaled to the 
current epoch (see Figure 21 in Fransson et al. 2013), 
and corresponds to complete ionization of the H and He 
mass within R s . 

The optical depth associated with the Hll region along 
the line of sight (los) can be estimated as 

T ff « 3.28 x 10~ 7 T^ 136 jy- 21 N^R S , (2) 

where T 4 = T e /(10 4 K), v is in GHz, N e is in cm" 3 and 
R s ~ J log dl is in pc. For T4 ~ 1, given the emission 
measure 4.6 x 10 6 < EM < 5.9 x 10 7 cm -6 pc, where 
EM = N^R S1 at frequencies 102 < v < 672 GHz the 



emission becomes nearly transparent as 2.28 x 10 < 
Tff < 1.2 x 10 -3 . The flux associated with the ionized 
component of the ejecta, can then be derived as 

2kT 

s " ff « r i\f 0)^-^ « * y-0 ' 1 - ( 3 ) 

Considering the solid angle Q subtended by the same 
radius, i? s , at all frequencies, for the lower limit Mhii ~ 
0.7 M©, i.e. N e = 6.8 x 10 3 cm" 3 , S u - (S„ - S Mm ) 
at 102 and 213 GHz, as Eq. 3 yields 0102 ff 00 3.9 mJy, 
#213// ~ 3.2 mJy, 5 3 45// ~ 3.1 mJy and S 672ff ~ 2.9 
mJy. These S Uff values (blue stars in in Figure 8) would 
well fit the SED (see blue fit in Figure 8). If the Hll 
density is considerably higher, as given by the upper limit 
Mhii ~ 2.5 M©, the derived S Uff fluxes would exceed the 
emission residuals by an order of magnitude. 

As a constraint to the free- free emission component, the 
hypothesized Hll region would also produce optical and 
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near-IR Hi recombination lines. The resultant Ha flux 
can be estimated as 



Sn a 



eN e N + ( -ttR) 



f 



And 2 ' 



(4) 



where e ex ola is the emissivity of the Ha line per 
unit volume, with a a = ola(T) = ^ n a n (T) the to- 
tal recombination coefficient, / is the volume filling fac- 
tor and d is the distance to the source. If one adopts 
e « 3.53 x 10- 25 T 4 " 0 92 erg cm" 3 s" 1 for T 4 ~ 0.5 
(e.g. Storey & Hummer 1995), the number of ions 
N+ = N e = 6.8 x 10 3 cm -3 as for the lower limit assumed 
for the density of the Hll region, and R s ~ 3.08 x 10 17 cm, 
for / ~ 1 one obtains Sn a ~ 1.26 x 10 -11 erg cm 2 s -1 . 
To match the Ha flux from the core measured by Frans- 
son et al. (2013) (see Figure 8 therein) at - 1.4 x 10" 14 
erg cm 2 s _1 on day 9000, the estimated Sn a would have 
to undergo 7.4 magnitudes of extinction. Such extinc- 
tion is possible but improbable. We note that the above 
Stf_ a also exceeds the Ha emission by the reverse shock, 
measured at ~ 2.0 x 10 -13 erg cm 2 s _1 on day 9000, 
and associated with a density 250 < N+ < 750 cm -3 
(Fransson et al. 2013). 

As regards the magnitude of the flux from the Br-7 line, 
the model estimate by Kjaer et al. (2010) on day 6840 
is a factor of ~ 10 -2 smaller than the Sn a derived from 
Eq. 4, while this is mainly associated with emission from 
the hot spots in the ER. 

In the absence of external X-ray heating, both the heat- 
ing and ionization would be powered by radioactive de- 
cays in the core. The resultant flux, estimated via models 
as in Kozma & Fransson (1998), would be several orders 
of magnitude smaller than the lower limit for S Vff . 

Another possible source of ionizing emission is a pulsar 
wind nebula (PWN) in the SNR interior. The properties 
of any PWN are very uncertain (see § 7), but there is the 
expectation that line emission would accompany free-free 
emission also in this case. 

5.2. Flat- Spectrum Synchrotron Emission 

A flat spectral index could also be attributed to a sec- 
ond synchrotron component. As shown in Figure 8, 
within 102 < v < 672 GHz, a synchrotron compo- 
nent with spectral index a p = —0.231 fits the residuals 
(S u — *S , 44 fit ) at 102 and 213 GHz. Synchrotron emis- 
sion with —0.4 < a < —0.1 fits the spectrum in the ra- 
dio/FIR transition and could originate from a compact 
source near the center of the SNR. In the case of a central 
pulsar, the synchrotron emission would be generated by 
the shocked magnetized particle wind (Gaensler & Slane 
2006). The scenario of a synchrotron-emitting PWN in 
the inner SNR is explored in § 7. 

5.3. Dust Emission 

If the excess emission is due to an additional synchrotron 
component, S u (a p ) oc ^~ 0,2 , this would provide a con- 
straint to the net dust emission, S Ud , in ALMA data. The 
subtraction S v —S±i m —S v (a p ), where S Vp = S 0p (u/iy 0 ) a p 
with So = 3.5 ± 3.1 mJy as from Table 2, leads to 
S 2 is d ~ P 0-7 m Jy, #345 d ~ 4.9 mJy and S 672d ~ 45.0 



mJy (black triangles in Figure 8). We note that #672 d 
coincides with the integrated flux density of the central 
feature of the related image, which extends over the inner 
ejecta as seen with HST (see Figure 4). The net dust can 
be fitted via a modified Planck curve of thermal emission, 



d 2 



k v B v (T)kvP B V (T), (5) 



where Md us t is the dust mass, k v = 3Q V /Apa v is the 
dust mass absorption coefficient, with Q v the absorption 
efficiency for spherical grains of density p and radius a Vl 
1 ^ P ^ 2 for insterstellar dust (e.g. Cortese et al. 
2014 and references therein) and B V (T) is the Planck 
function. As shown in Figure 8, the best fit of both 
ALMA fluxes, as reported in this paper, and Herschel 
fluxes, as from Matsuura et al. (2011), yields f3 = 1.23 
and dust temperature T — 22.2 K. The thermal peak 
of the SED has been previously fitted (Lakicevic et al. 
2011, 2012a, b; Matsuura et al. 2011) with temperatures 
estimated between 17 and 26 K, while Indebetouw et al. 
(2014) fit amorphous carbon dust at T = 26 ±3 K. 

The sum of the main synchrotron component, 
S y oc ^ _0 - 73 , and the emission component from 
dust grains at T ~ 22 K is lower than the measured 
ALMA flux densities at 213 and 345 GHz (see dashed 
gray fit in Figure 8). To match the emission excess 
of ~ 3 mJy in this frequency range, we could also 
postulate a second dust component. This would require 
very cold dust at T < 5 K, i.e. at temperatures where 
the assumption of either amorphous carbon or silicates 
leads to dust masses implausibly large for physically 
realistic grains. In particular, as from Eq. 5, at 345 
GHz a flux density of ~ 3 mJy would require dust at 
T ~ 3 K with Mdust ~ 50 M 0 , as obtained by using 
k v ~ 2.5 cm 2 g _1 for amorphous carbon (Zubko et al. 

1996, 2004). We note that warmer dust at T - 180 
K has been identified by Dwek E. et al. (2010) in the 
ER, where the dust grains are likely collisionally heated 
by the expanding radiative shocks (Bouchet et al. 2006). 

6. SPECTRAL INDEX VARIATIONS 

The spectral index distribution across the remnant is in- 
vestigated via multi-frequency spectral index maps and 
the T-T plot method (Costain 1960; Turtle et al. 
1962). 

Spectral index maps are derived from images at 44, 94, 
102, 213 and 345 GHz. In Figure 9, the maps resulting 
from the combination of images at two, three, four and 
five frequencies are shown. All images are reduced with 
identical procedure and restored with a circular beam 
of Of! 8. ATCA data at 44 and 94 GHz are scaled to 
day 9280, via exponential fitting parameters derived for 
ATCA flux densities from day 8000, as measured at 8.6 
and 9 GHz (Zanardo et al., in preparation). While the 
two-frequency maps are derived from direct division of 
the flux densities, the maps derived from 3, 4 and 5 fre- 
quencies are more accurate, as the solution minimizes the 
function 
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Figure 9. In color scale are maps of the spectral index distribution in SNR 1987A, as derived from images at 44, 94, 102, 213, and 345 
GHz. The source images are reduced with identical procedure and restored with a circular beam of Of/ 8. To match the average epoch 
of ALMA observations, ATCA data at 44 and 94 GHz are scaled to day 9280, via exponential fitting parameters derived for ATCA flux 
densities from day 8000, as measured at 8.6 and 9 GHz (Zanardo et al., in preparation). The spectral maps are derived from flux densities 
at 2, 3, 4, and 5 frequencies, and labelled as 2f, 3f, 4f-> an d 5/, respectively (top right corner). The frequencies used in each map are 
indicated, in GHz, on the top left corner. Image regions below the highest rms noise level are masked. The upper limit of the spectral 
index color scale is set to a = 0.25, therefore map regions with spectral indices greater than 0.25 appear in white. In grey scale are maps 
of the error associated with a in the power-law fit, S u ~ used for spectral maps of flux densities at 3, 4, and 5 different frequencies, 
and labelled as 3f, 4f> an< l 5f, respectively. To emphasize the error distribution, while using a linear gradient of the grey color scale, map 
regions with errors greater that 0.01 appear in black. In the top right map, the contours of the 44 GHz image, resolved with a 0 / /25 circular 
beam, are shown at the 15% and 60% emission levels (in blue). The red cross indicates the VLBI position of SN 1987A as determined by 
Reynolds et al. (1995) [RA 05 h 35 m 27!968, Dec -69° 16' ll 7 / 09 (J2000)]. 
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Figure 10. Top — Four-frequency (44-94-213-345 GHz) spectral index intensity map with an angular resolution of 0'/8, where ATCA 
data at 44 and 94 GHz are scaled to day 9280. The map is superimposed with the contours of the 44 GHz image, resolved with a Of/ 25 
circular beam, at the 15% and 60% emission levels (in blue). The regions used to derive the spectral indices oltt-, via T-T plots, are 
indicated by grey/white squares. Bottom — The T-T plots are applied to 6 pairs of frequencies: 44 — 102 GHz, 102 — 213 GHz, 213 — 345 
GHz, 44 — 94 GHz, 44 — 213 GHz, and 44 — 345 GHz. The controlling beam associated with each pair of frequencies is shown on the lower 
left corner of the images. The T-T boxes are color-coded according to the resultant oltt value. For image pairs of angular resolution lower 
than 0'/8, larger T-T boxes are used as outlined. All T-T boxes are drawn on top of the 4-frequency spectral map, highlighted by contours 
spaced at Aa = 0.1, which are colored according to the associated spectral index. The red cross indicates the VLBI position of SN 1987A 
as determined by Reynolds et al. (1995). The oltt values resulting from the six frequency pairs are detailed in Table 4. 
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Table 4 

Regional Spectral Indices 











fo a \( a ) ( e ) 






OLM^ l^l,2,3,4j 


Region (c) 


^. (d) 

Fit 


(44,102) 


(102,213) 


(213,345) (44,94) 


(44,213) 


(44,345) 


(44,94,213,345) 



Nl 



N2 



y 

ODR 

Abs 
(0,0) 
V 

ODR 

Abs 
(0,0) 
V 

ODR 
Abs 
(0,0) 



-0.33 ±0.02 
-0.41 ±0.02 
-0.40 ±0.02 
-0.05 ±0.04 



-0.38 ±0.03 



-0.37 ±0.03 
-0.43 ±0.05 



-0.28 ±0.08 
-0.31 ±0.09 
-0.27 ±0.08 
-0.36 ±0.14 
-0.73 ±0.03 
-0.76 ±0.03 
-0.77 ±0.03 
-0.84 ±0.04 



-0.39 ±0.02 
-0.42 ±0.02 
-0.35 ±0.03 
-0.40 ±0.04 
-0.81 ±0.08 



-0.89 ±0.07 



-0.39 ±0.09 



-0.85 ±0.14 



E0 



El 



E2 



y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 



-0.95 ±0.02 
-0.95 ±0.01 
-0.95 ±0.02 
-0.97 ±0.02 



-0.07 ±0.02 
-0.07 ±0.03 
-0.06 ±0.03 
-0.03 ±0.03 



-1.05 ±0.04 

-1.12 ±0.04 
0.11 ±0.07 



-0.48 ±0.02 
-0.60 ±0.05 
-0.47 ±0.02 
-0.46 ±0.02 



-0.55 ±0.03 


-0.68 ±0.03 


-0.60 ±0.03 


-0.59 ±0.03 


-0.72 ±0.04 




-0.54 ±0.03 


-0.66 ±0.03 




-0.69 ±0.07 


-0.71 ±0.06 




-0.65 ±0.04 


-0.66 ±0.04 


-0.52 ±0.13 


-0.68 ±0.04 


-0.68 ±0.04 




-0.70 ±0.05 


-0.65 ±0.04 




-0.78 ±0.08 


-0.68 ±0.04 




-0.60 ±0.03 


-0.73 ±0.03 


-0.55 ±0.12 


-0.63 ±0.03 


-0.74 ±0.02 




-0.65 ±0.03 


-0.72 ±0.03 




-0.74 ±0.04 


-0.73 ±0.03 





SI 



y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 



-0.05 ±0.05 
-0.10 ±0.04 
-0.23 ±0.03 
-0.22 ±0.03 



-0.45 ±0.02 

-0.44 ±0.02 
-0.53 ±0.07 



-0.61 ±0.02 

-0.60 ±0.02 
-0.95 ±0.18 



-0.52 ±0.07 

-0.56 ±0.07 
-0.67 ±0.09 



-0.63 ±0.07 



W 



WO 



Wl 



W2 



y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 



-0.77 ±0.02 
-0.78 ±0.02 
-0.77 ±0.02 
-0.76 ±0.02 



-0.11 ±0.03 
-0.11 ±0.02 
-0.10 ±0.03 
-0.18 ±0.03 



0.25 ± 0.04 
0.17 ±0.02 
0.08 ± 0.04 
0.16 ±0.03 



-0.40 ±0.06 
-0.43 ±0.06 
-0.43 ±0.06 
-0.49 ±0.05 



-0.69 ±0.03 


-0.49 ±0.09 


-0.53 ±0.05 


-0.77 ±0.04 


-0.57 ±0.10 






-0.54 ±0.10 




-0.66 ±0.03 


-0.49 ±0.11 




-0.93 ±0.09 


-0.61 ±0.09 


-0.37 ±0.11 


-0.97 ±0.08 






-0.99 ±0.08 


-0.60 ±0.09 




-0.75 ±0.10 


-0.51 ±0.15 




-0.40 ±0.11 


-0.50 ±0.02 


-0.47 ±0.11 


-0.44 ±0.11 


-0.52 ±0.02 




-0.35 ±0.12 


-0.49 ± 0.02 




-0.52 ±0.14 


-0.51 ±0.02 





CI 



y 

ODR 

Abs 
(0,0) 

y 

ODR 

Abs 
(0,0) 



-0.86 ±0.08 
-0.90 ±0.08 
-0.86 ±0.08 
-0.84 ±0.08 



-0.00 ±0.04 
-0.00 ±0.03 
-0.00 ±0.03 
-0.09 ±0.04 



1.20 ±0.04 
1.97 ±0.09 
2.25 ±0.11 



-0.43 ±0.02 

-0.44 ±0.01 
-0.43 ±0.02 



-0.40 ±0.08 

-0.47 ±0.08 
-0.73 ±0.15 



-0.01 ±0.06 



a The spectral indices oltt are derived from 5 images at 44, 94, 102, 213, and 345 GHz. The images are analyzed in 6 frequency pairs, 
as indicated in Figure 10. 

b The median spectral index, olm^ •> is derived by fitting a Gaussian to the histogram of the olm values resulting in the spectral map 
obtained from images at 44, 94, 213, and 345 GHz, in the corresponding T-T region (see top map in Figure 10). 
c The regions selected for the T-T plots are designated in Figure 10. 

d The given spectral indices are derived from different linear interpolations: least-square fit from vertical squared errors, a y ; orthogonal 
distance regression, qlodr', robust fitting from the square root of absolute residuals, a Abs', linear regression with forced zero interception 
(5 , i=0,5 2 =0),a ( o j o)- 

e The errors on oltt and olm are the 1 a error on the slope of the linear fit combined with the uncertainty in the flux calibration. 
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Sj (mjy/beam) 



Figure 11. Temperature-temperature (T-T) plot of the whole 
SNR, where flux densities (Si vs S2) are plotted instead of 
brightness temperatures. Four frequency pairs are included: 
44-94 GHz (green hollow squares); 44-102 GHz (solid black 
squares); 102-213 GHz (solid blue squares); 213-314 GHz (solid 
red squares). Data points for the 44-94 GHz frequency pair are 
not fitted. The given spectral indices are derived from different 
linear interpolations: least-square fit from vertical squared errors, 
cty (blue); orthogonal distance regression, olodr (green); robust 
fitting from the square root of absolute residuals, ctAbs (magenta); 
linear regression with forced zero interception (Si = 0,^2 = 0), 
«(o,o) ( red )- 

at each pixel of coordinates (p, q) : for images at fre- 
quencies, with n > 3. 

In Figure 9, it can be seen that the spectral indices over- 
all become flatter as the frequencies reach the FIR. For 
all maps, the spectral index, ckm? varies between —1.25 
and 0.25, while in the 5— frequency map the index range 
narrows to —0.96 < olm < 0.18. In most of the multi- 
frequency maps, olm is steeper on the eastern half of 
the SNR, while flat spectral indices surround the cen- 
ter, where the bulk of dust sits (see § 3), and extend 
onto the NW and SW regions of the remnant. With 
reference to the 5— frequency map, — 0.9 < olm < — 0.6 
on the eastern lobe and —0.4 < < 0 on the west- 
ern side of the SNR, with olm ~ 0 around the central 
region and —0.4 < olm ^ —0.3 predominantly in the 
NW quadrant (PA - 330°) and in the SW quadrant (PA 
- 210°). 

The T-T plot method is applied to two images at dif- 
ferent frequencies, each obtained with identical reduc- 
tion process and with the same angular resolution. The 
spectral variations are assessed over image regions not 
smaller than the beam size, where the spectral index, 
oltt j is determined from the flux density slope m, where 
S2 = m S\ + q. The regions used for the T-T plots are 
shown in Figure 10. Six frequency pairs are considered: 
44-102 GHz, 102-213 GHz, 213-345 GHz; and 44-94 



GHz, 44-213 GHz and 44-345 GHz. Different box sizes 
are used to suit the controlling beam of each frequency 
pair. All derived 

OLTT 

values are listed in Table 4. 

Similarly to the trend of oltt values become flatter 
at higher frequencies (see Figure 11). From 102 to 213 
GHz, oltt ~ —0.1 across the whole remnant (see Fig- 
ure 10), while from 213 to 345 GHz the spectral distribu- 
tion appears markedly split in two larger regions (see also 
the two- frequency map in Figure 9), with very flat indices 
on the western half of the SNR and steep indices on the 
eastern side. For the frequency pairs 102-213 and 213- 
345 GHz, the T-T plots for the eastern region (see E in 
Table 4) give Aol T t = ^(102, 213) -«(2i3, 345) = 1.0±0.07. 
This could be indication of a local spectral break at 
213 GHz. Using v c = 213 GHz in Eq. 1, for B = 20 
mG, r e ~ 35 yr. However, given the high Mach num- 
ber (A4 ~ 10 4 ) of the eastbound shocks (Zanardo et 
al. 2013), likely sub-diffusive particle transport (Kirk 
et al. 1996) by the shock front and, consequently, local 
magnetic-field amplifications (Bell 2004), it is possible 
that the CR in the eastern lobe, radiating at v c , are al- 
ready past their synchrotron lifetime. 

Flat spectral indices in the western lobe extend north 
and south at both 213-345 GHz and 44-213 GHz (see 
Figure 10), while T-T plots from higher resolution images 
show a narrower north-south alignment of the flat regions 
(see Nl, CI and W2 in Figure 10 and Table 4), with 
-0.5 < ol tt < -0.3. The T- T plots for 44-94 GHz and 
44-345 GHz also yield the steepest spectral indices in 
region Wl, this might be due to the local emission drop 
in the NW sector of the ER, visible in the images at 94 
and 345 GHz at PA - 300°. 

As discussed in § 5, the flat-spectrum western regions 
could be linked to a PWN. We note that spectral maps 
of the Crab Nebula via observations centered at 150 GHz 
(Arendt et al. 2011) have identified spectral indices 
around —0.2 in the inner central regions of the PWN, 
while spectral indices ~ —0.3 have been associated with 
the PWN periphery. 

7. PWN CONSTRAINTS 

As discussed in § 3.1, the emission at 102 and 213 GHz 
in the I44— subtracted images appears to peak west of 
the SN site (Reynolds et al. 1995) (Figures 1, 2), and 
to mainly extend west of the optical ejecta (see Fig- 
ure 5). Besides, both the spectral maps and T-T plots 
(Figures 9, 10) show that flat spectral indices can be asso- 
ciated with the center- west regions of the SNR (see § 6). 
These results could be explained by a possible PWN, 
powered by a pulsar likely located at a westward offset 
from the SN position. 

The pulsar-kick mechanism has been linked to asymme- 
tries in the core collapse or in the subsequent supernova 
explosion, presumably due to asymmetric mass ejection 
and/or asymmetric neutrino emission (Nordhaus et al. 

2012; Podsiadlowski 2005; Wongwathanarat et al. 
2013). As evidence for the natal kick, neutron star (NS) 
mean three-dimensional (3D) speeds have been estimated 
at v NS ~ 400 ± 40 km s _1 (Hobbs et al. 2005), while a 
transverse velocity of ~ 1083 km s _1 has been detected 
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by Chatterjee et al. (2005), which would imply a 3D 
NS birth velocity as high as 1120 km s _1 (Chatterjee et 
al. 2005). In the context of SNR 1987A, by day 9280 
the NS could have travelled westwards of the SN site by 
~ 20 — 80 mas, while for an impulsive kick of the same 
order as ^ns a distance of ~ 42 ± 5 mas would have been 
covered. With a western offset of ~ Of! 05 , the NS would 
be situated inside the beam of the subtracted im- 
ages associated with the peak flux density, both at 102 
and 213 GHz. If we take into account the ~60 mas un- 
certainty intrinsic to image alignment (see § 2), as well 
as the error of 30 mas in each coordinate of the SN VLBI 
position (Reynolds et al. 1995), the NS could be located 
near the emission peak as seen in the subtracted im- 
ages at 102 and 213 GHz. We note that the inner fea- 
ture of fainter emission detected in the SNR at 44 GHz, 
as aligned with VLBI observations of the remnant (Za- 
nardo et al. 2013), is centered ~ 60 mas west of the SN 
site. 

If a pulsar is embedded in the unshocked eject a, the 
PWN would be in its early stages of evolution, likely sur- 
rounded by uniformly expanding gas (Chevalier & Frans- 
son 1992). Diffuse synchrotron emission from the PWN 
would be due to the relativistic particles, produced by 
the pulsar, accelerated at the wind termination shock 
(Kirketal. 2009). 

Assuming a power-law energy distribution of electrons, 
i.e. the particle density N e is expressed as N e (E) oc 
KE~ S , where s = 1 — 2a and K oc (m e c 2 ) s_1 , the syn- 
chrotron emission of the PWN can be written as 



field, can be simplified as 



KB 



!(* + !) 1(1 -s) 



(7) 



where B rwN is the nebular magnetic field strength. Not- 
ing that, in radio observations, the energy in electrons 
cannot be separated from that in the magnetic field 
(Reynolds et al. 2012), the equipartition magnetic-field 
strength could be derived as (e.g. Longair 2011; see re- 
vised formula by Arburtina et al. 2012; Beck & Krause 
2005) 



G 0 G(/C + 1) 



(1-8) 



fdO* 



(5+5) 



(8) 



where Go is a constant, G = G(y, s) is the product of 
different functions varying with the minimum and max- 
imum frequencies associated with the spectral compo- 
nent and the synchrotron spectral index (Beck & Krause 
2005; Longair 2011), K is the ion/electron energy ra- 
tio, / is the volume filling factor of radio emission, and 
# PWN = R pwN /d is the angular radius. Considering 
0"05 < R PWN < 0715, -0.4 < a < -0.1 (as from § 5), 
102 < v < 672 GHz, S v « 3 mJy, and taking K « 100 
(Beck & Krause 2005) while / « 0.5, Eq. 8 leads to 
1 < ^pwn < 7 mG. For 0.2 </ < 0.9, 2 < £ PWN < 5 
mG. Since the equipartition is a conjecture for young 
SNRs and no longer valid when the spectral index is 
flatter than —0.5, these B PWN estimates might be inac- 
curate. 

The energy inside the PWN, due to the PWN magnetic 



} W~^pwnW : 



Bl 



8tt 



(9) 



where the magnetic field is considered uniform 
and isotropic inside the PWN volume, V PWN (t) = 
4/3ffl 3 PWN (i). For 0"05 < i? PWN < C/15, and, as for 
the parameters used in Eq. 8, 1 < ^ PW n ^ ? m G, 
at t = £25 = 8.0 x 10 8 s (« 9280 days), we estimate 
0.9 x 10 43 < £ PWN;B < 1.2 x 10 46 erg. According to 
models by Chevalier & Fransson (1992), about 30% of 
the total energy input into the PWN, E PWN , goes to the 
internal magnetic pressure in the PWN, while most of 
the remaining pulsar spin-down energy would drive the 
PWN expansion into the ejecta. 

In terms of integrated radio luminosity, calculated 

as 



^rad(^) = 



(10) 



if z^min — 102 GHz and ^ ma x — 

672 GHz bracket the 
frequency range in which the PWN is detected, S u ~ 3 
mJy leads to L rad « 5.4 x 10 33 erg s _1 . The derived 
L ra d is comparable with the limit of L opt < 5 x 10 33 erg 
s _1 given by Graves et al. (2005) for a compact source 
in the optical band from 290 to 965 nm at t = 6110 
days. A similar limit has been placed on the 2-10 keV 
X-ray luminosity, Lx < 5.5 x 10 33 erg s _1 , using Chandra 
images (Park et al. 2004). Since these L v estimates are 
upper limits and free-free emission may be a substantial 
component of the radio luminosity (see § 5.1), we can 
take L ra d « 10 33 erg s _1 as a realistic upper limit. 

In the free-expansion regime (Chevalier & Fransson 
1992), rsj 1% of the pulsar power is emitted by the 
shock wave in the supernova and additional emission 
from the pulsar nebula is expected. Given this, we can 
set i^spin ~ 10 35 erg s _1 as an upper limit for the spin- 
down power of the pulsar. For a typical pulsar sur- 
face dipole magnetic field B s ~ 10 12 G (Manchester 
et al. 2005b), this spin-down luminosity corresponds 
to a pulsar period P ~ 0.15 s and characteristic age 
r c = Pj (2P) rsj 10 5 years. Lower luminosities would 
imply lower dipole magnetic fields and/or longer pulsar 
periods. Such parameter ranges are plausible for the pu- 
tative pulsar at the centre of SN 1987A, as there is good 
evidence that many pulsars are born with a spin period 
not much different to their present period (Gotthelf et 
al. 2013; Popov & Turolla 2012). 

Following Chevalier (1977), the velocity at the outer 
edge of the nebula can be defined as 



= 6/125 Es P int 
5V132 irpejt 3 



(11) 



assuming the PWN is freely expanding and has constant 
density p e j. At t = £25, setting 10 -18 < 



< 10 



-19 



cm -3 in the central region of the SNR, as from the den- 
sity model by Blinnikov et al. (2000) (see Figure 21 in 
Fransson et al. 2013), the swept-up shell velocity be- 
comes 260 < v rwN < 410 km s _1 , which leads to R PWN 
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Figure 12. Dispersion measure (DM) associated with a 
possible pulsar in SNR 1987A, assuming an ionized fraction 
of the ejecta 0.1 Mq < Mhii < 2.5 Mq, uniformly distributed 
within a spherical region with radius 0.01 < i^nil S 0-10 pc (i.e. 
CK'05 < Run i$ 0 / /42). The color scheme changes from white to 
black for increasing DM values, which range from 100 to 6000 
cm -3 pc, as indicated by the contour labels. The contour at DM 
=1000 cm -3 pc is highlighted in magenta. 

not greater than 0'/05. Since this is well below the resolu- 
tion of the ATCA and ALMA images presented here, the 
emission from a possible PWN would appear as a point 
source. 

As mentioned in § 5, a pulsar embedded in the SNR 
interior would emit ionizing radiation within the inner 
layers of the ejecta. While an X-ray pulsar has yet to 
be detected (Helder et al. 2013), illumination of the 
inner ejecta by X-ray flux has been reported by Larsson 
et al. (2011) though attributed to the reverse and re- 
flected shocks, as well as to shocks propagating into the 
ER. 

Given the early stages of a possible PWN, the 
ionized ejecta would be mainly due to the radi- 
ation from the various shocks. If the Hll region 
is assumed to be spherical (see § 5), the re- 
lated dispersion measure (DM) can be derived as 
DM(M H ii,#hii) « 3M H ii/(m p 47ri^ n ). Neglecting 
clumping in the ejecta, for 0.1 M 0 < Mhii ^ 2.5 M 0 
and 0.01 < R mi < 0.10 pc (i.e. Of! 05 < Rmi < 0 / /42), 
the resulting DM is shown in Figure 12. 

8. CONCLUSIONS 

We have presented a comprehensive morphological and 
spectral analysis of SNR 1987A based on both ATCA and 
ALMA data (Cycle 0), from 1.4 to 672 GHz. We have 
investigated the components of the SNR emission across 
the transition from radio to FIR, as the combination of 
non-thermal and thermal emission. A summary of our 



findings is as follows: 

1. To decouple the non-thermal emission from that 
originating from dust, the synchrotron component, 
as resolved with ATCA at 44 GHz, and the dust 
component, as imaged with ALMA at 672 GHz, 
were subtracted from the datasets at 94, 102, 213, 
345 and 672 GHz. The images derived from the 
subtraction of the scaled model flux density at 672 
GHz, highlight the ring-like synchrotron emission 
that originates over the ER. The images similarly 
obtained from the subtraction of the 44 GHz model, 
show residual emission mainly localized west of the 
SN site. 

2. Analysis of the emission distribution over the ER, 
in images from 44 to 345 GHz, shows a gradual 
decrease of the east-to-west asymmetry ratio with 
frequency. We attribute this to the shorter syn- 
chrotron lifetime at high frequencies. We estimate 
that, at frequencies higher than 213 GHz, the elec- 
trons might be unable to cross the eastern emission 
region. 

3. Across the transition from radio to FIR, the SED 
suggests additional emission components beside the 
synchrotron main component (S u oc u~ 0 - 73 ) and 
the thermal emission originating from dust grains 
at T ~ 22 K. We argue that this excess emission 
could be due to a second flat-spectrum synchrotron 
component with —0.4 < a < —0.1. This could 
imply the presence of a PWN originating from an 
embedded pulsar. 

4. Spectral index measurements across the SNR, from 
44 to 345 GHz, show predominantly flat spectral 
indices, —0.4 < a < —0.1, in the western half of the 
remnant at frequencies above 102 GHz, while a ~ 0 
around the central region. From 102 to 213 GHz, 
— 0.1 < a < 0 across the whole remnant. From 213 
to 345 GHz, the steepening of the spectral indices 
over the eastern lobe might be indication of a local 
spectral break. 

5. Results from both the morphological and spectral 
analysis might suggest the presence of a PWN in 
the SNR interior, powered by a pulsar likely located 
at a westward offset from the SN position. In this 
scenario, a NS embedded in the unshocked ejecta 
might be emitting diffuse synchrotron emission 
S vp ~ 3 mJy above 102 GHz. The possible PWN 
would have magnetic field strength 1 < B PWN < 7 
mG. The upper limit of the integrated radio lumi- 
nosity is derived as L ra d ~ 5.4 x 10 33 erg s _1 , which 
leads to an upper limit of the pulsar spin-down 
power £? S pin ~ 10 35 erg s -1 . While effects such 
as luminosity, beaming, scattering and absorption 
might have prevented the detection of pulsed emis- 
sion in the SNR for over 20 years (Manchester 

1988, 2007), a renewed pulsar search with the 
Parkes telescope is currently in progress. Future 
observations with ALMA, ATCA and VLBI are 
needed to probe the SNR for a possible PWN. 
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